Criticality of the "critical state" of granular media: 
Dilatancy angle in the Tetris model 
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The dilatancy angle describes the propensity of a granular medium to dilate under an applied shear. 
Using a simple spin model (the "Tetris " model) which accounts for geometrical "frustration" 
effects, we study such a dilatancy angle as a function of density. An exact mapping can be drawn 
with a directed percolation process which proves that there exists a critical density pc above which 
the system expands and below which it contracts under shear. When applied to packings constructed 
by a random deposition under gravity, the dilatancy angle is shown to be strongly anisotropic, and 
it constitutes an efficient tool to characterize the texture of the medium. 



I. INTRODUCTION 

Granular materials p| give rise to a number of original 
phenomena, which mostly result from their peculiar rhe- 
ological behavior. Even using the most simple descrip- 
tion of the grains (rigid equal-sized spherical particles) 
a granular system displays a rather complex behavior 
which shows that the origin of this rheology has to be 
found at the level of the geometrical arrangement of the 
grains. 

Guided by these considerations, models have been pro- 
posed to account for the geometrical constraints of as- 
semblies of hard-core particles |§-^. The motivation 
of these models is not to reproduce faithfully the local 
details of granular media, but rather to show that sim- 
ple geometrical constraints can reproduce under coarse- 
graining some features observed in real granular media. 
Along these lines, one of the most impressive examples 
is the "Tetris" model [g,|| which, in its simplest version, 
is basically a spin model with only hard core repulsion 
interactions. This model has been introduced in order to 
discuss the slow kinetics of the compaction of granular 
media under vibrations. In spite of the simplicity of the 
definition of the model, the kinetics of compaction has 
been shown to display a very close resemblance to most 
of the experimentally observed features of compaction 0] 
and segregation pO| ]. 

Our aim is here to consider again the Tetris model 
and to focus on a basic property of the quasi-static 
shearing of a granular assembly. It is well known since 
Reynolds that dense granular media have to dilate in or- 
der to accommodate a shear [Q, whereas loose systems 
contract. This observation is important since it gives ac- 
cess to one of the basic ingredients (the direction of the 
plastic strain rate) necessary to describe the mechanical 



behavior in continuum modeling. The dilatancy angle 
is defined as the ratio of the rate of volume increase to 
the rate of shearing. Denoting with e^y the component 
xy of the strain tensor e, Fig.(|l|) illustrates an experi- 
ment where a shear rate ixy is imposed together with a 
zero longitudinal strain rate ixx = 0, and the volumet- 
ric strain rate (here vertical expansion) iyy is measured. 
The direction of the velocity of the upper wall makes an 
angle ip with respect to the horizontal direction. This 
angle is called the dilatancy angle, ^. In this particular 
geometry we have 



tan(^) = ^ 



(1) 



More generally, the tangent of the dilation angle is the 
ratio of the volumetric strain rate (tr(£)) to the deviatoric 
part of the strain rate. 




FIG. 1. 

Numerous experimental studies have confirmed the va- 
lidity of such a behavior, and have lead to extensions such 
as what is known in soil mechanics as the "critical state" 



concept iQ. Assuming that the incremental (tangent) 
mechanical behavior can be parametrized using only the 
density of the medium, p, a loose medium will tend un- 
der continuous shear towards a state such that no more 
contraction takes place, i.e. it will assume asymptoti- 
cally a density pc such that '(/'(Pc) — 0. This state is by 
definition the "critical state" . Conversely, if the strain 
were homogeneous, a dense granular media would dilate 
until it reached the critical state density pc- However, 
for dense media, the strain may be localized in a nar- 
row shear band which may allow a further shearing with- 
out any more volume change so that the mean density 
may remain at a value somewhat higher than the critical 
value. Recent triaxial tests |1^ in a scanner apparatus 
have however shown that in the shear bands the density 
of the medium was quite comparable to the critical den- 
sity, thus providing further evidence for the validity of 
the critical state concept. 

The word "critical" used in this context has become 
the classical terminology, but it has no a priori relation 
with any kind of critical phenomenon in the statistical 
physics vocabulary [Q . One of the results presented in 
this article is to show that indeed the critical state of soil 
mechanics is also a critical point in the sense of phase 
transitions, for the Tetris model considered here. 



II. MODEL AND DEFINITION OF DILATANCY 

A group of lattice gas models in which the main ingre- 
dient is the geometrical frustration has been introduced 
recently under the name Tetris P,pi. 




FIG. 2. 

The Tetris model is a simple lattice model in which 
the sites of a square lattice can be occupied by (in its 
simplest version) a single type of rectangular shaped par- 
ticle with only two possible orientations along the princi- 
pal axis of the underlying lattice. A hard core repulsion 
between particles is considered so that two particles can- 
not overlap. This forbids in particular that two nearest 



neighbor sites could be both occupied by particles aligned 
with the inter-site vector. An illustration of a typical ad- 
missible configuration is shown schematically in Fig.(@). 
More generally one can consider particles that move on 
a lattice and present randomly chosen shapes and sizes 
Q. The interactions in the system obey to the general 
rule that one cannot have particle overlaps. The inter- 
actions are not spatially quenched but are determined 
in a self-consistent way by the local arrangements of the 
particles. 

The definition of the dilation angle as sketched in 
Fig.(n^) is difficult to implement in practice in the 
Tetris model due to the underlying lattice structure 
which defines the geometric constraints only for particles 
on the lattice sites, and not in the continuum. 

We may however circumvent this difficulty through the 
following construction illustrated in Fig. (|^) . Wc consider 
a semi-infinite line starting at the origin and oriented 
along one of the four cardinal directions. This line is 
(and all the sites attached to it are) pushed in one of the 
principal directions of the underlying square lattice by 
one lattice constant. In the following, we will consider 
only a displacement perpendicular to the line, although 
a parallel displacement may also be considered. As this 
set of particles is moved, all other particles which may 
overlap with them are also translated with the same dis- 
placement. In this way, we determine the set T) of par- 
ticles which moves. In the sequel, we will show that this 
domain is nothing but a directed percolation cluster p^ ] 
grown from the line. Anticipating on the following, the 
mean shape of T) will be shown to be a wedge limited by a 
generally rough boundary whose mean orientation forms 
an angle -0 with the direction of motion. The angle tp 
can be shown to be exactly equal to the dilatancy angle 
as defined previously. 




FIG. 3. 

Exploiting the non-overlap constraint, we may simply 
determine the rule for constructing the domain D. Let us 
choose the particular case of a displacement in the direc- 
tion (1,0), and consider a non-empty site {i,j) which is 



displaced. The particles which may have to be displaced 
together with site {i,j) can be identified easily: 
-If the particle in (i,j) is horizontal: 

• (i + 1, j) if the site is occupied by a particle with any 
orientation. 

• (« + 2, j) if the site is occupied by a horizontal particle. 
- If the particle in (i,j) is vertical: 

• (« + 1, j) if the site is occupied by a particle with any 
orientation. 

• (i + 1, j ± 1) if the site is occupied by a vertical particle. 
Using these rules, it is straightforward to identify the 

cluster of particles V. The model thus appears to be a 
directed percolation problem with a mixed site/bond lo- 
cal formulation. Thus unless long range correlations are 
induced by the construction of the packing, the resulting 
problem will belong to the universality class of directed 
percolation. The density of particles, p € [0,1], in the 
lattice plays the role of the site or bond presence proba- 
bility, i.e. the control parameter of the transition. 

Let us recall, for sake of clarity, some properties of 
the two-dimensional directed percolation. For p < pc 
(where pc is the directed percolation threshold) , a typical 
connected cluster extends over a distance of the order of 
^/f in the parallel direction (the preferential direction) 
and a distance ^_l in the perpendicular direction. For 
p > Pc there appears a directed percolating cluster which 
extends over the whole system. This cluster possesses a 
network of nodes and compartments. Each compartment 
has an anisotropic shape similar to the connected cluster 
below Pc, characterized by ^// in the parallel direction 
and f j_ in the perpendicular direction. On both sides 
of the percolation transition, the two lengths present the 
power-law behavior ^// ^ \p—Pc\~'^^^ and^j_ ^ \p—pc\~'^-^- 

III. MONOCRYSTAL 

Let us first examine a simple geometrical packing. 
There exist (two) special ordered configurations of par- 
ticles such that the density can reach unity (one particle 
per site). This corresponds to a perfect staggered dis- 
tribution of particle orientations. Thus a simple way of 
continuously tuning the density is to randomly dilute one 
of these perfectly ordered states. In this case, if a site 
is occupied by a particle, its orientation is prescribed. 
Therefore the above rules can be easily reformulated as 
a simple directed site percolation problem in a lattice 
having a particular distribution of bonds (up to second 
neighbors). Fig.(0) illustrates the specific distribution of 
bonds corresponding to such an ordered state. 




For p = 1, suppose that the initial seed is (0,j) for 
j > and this line is pushed in the x direction. Then 
the infinite cluster is the set of sites {i,j) such that j > 
—i, for a vertical spin at the origin. Thus moving the 
semi-infinite line (seed) introduces vacancies in the lattice 
which was initially fully occupied. The system dilates 
and its dilation angle is ipi = 7r/4. 

As p is reduced, the orientation of the boundary 
changes up to the stage where it becomes parallel to the 
X axis for p = p^. At this point the dilatancy is zero. 
A motion is possible without changing the volume. This 
point corresponds precisely to the directed percolation 
threshold (using the precise rules defined above). 

From the theory of directed percolation, we can di- 
rectly conclude that the behavior of the dilatancy angle 
Tp in the vicinity of pfi obeys 



tan('(/)) (X {p~ pr) 



7/-^i 



(2) 



where the correlation length exponents are v/j « 1.732 
and vx^ « 1.096 independently of the lattice used. 

A further decrease of p leads to a subcritical regime 
where only a finite cluster is connected to the initial seed. 
Only a finite layer of thickness ^^ ex {pn — p)~'^'' along 
the y-axis is mobilized. This means that it not possible 
to define in the same way the dilation angle for p < pR 
(negative angles). What happens in practice is that for 
p < PR the shearing produces a compaction of the sys- 
tem in front of the semi-infinite line pushing the system. 
Fig. (H) summarizes schematically the situation for all the 
values of p. The horizontal line corresponds to p = pr 
and a zero dilation angle. 
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We performed numerical simulations of this problem 
using a transfer matrix algorithm which allowed to gen- 
erate system of size up to 10'* x 3- 10*. These large system 
sizes allowed for a very accurate determination of the di- 
latancy angle as a function of the occupation probability 
(density) p. Fig. (||) shows the boundaries of the domains 
Dior p — 0.58, close to the directed percolation threshold 
PR, and p = 0.7. 






(a) 



(b) 



FIG. 6. 



Fig.(R) shows the estimated dilatancy angle as a func- 
tion of the density of particles. Angles are evaluated on 
lattice of size 10* x 3- 10"* and are averaged over 100 differ- 
ent realizations. The onset of dilatancy is thus estimated 
to be 



Pi^ = 0.583 ± 0.001 



(3) 
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The singular variation of ip close to the onset of dila- 
tancy Eq. has been checked to be consistent with our 
numerically determined values as shown by the dotted 
curve in Fig.([7|) which corresponds to the expected criti- 
cal behavior. 



IV. RANDOM SEQUENTIAL DEPOSITION 

It is worth emphasizing that the directed percolation 
problem associated with the dilatancy angle determina- 
tion is simply a site percolation problem in the above 
special case where each site is assigned only one possible 
orientation for the particle. In the more general case, the 
way the cluster is grown locally depends on the specific 
orientation of the particle. Thus it is a mixed site/bond 
percolation problem. Therefore, depending on the way 
the system has been built, the onset for dilatancy, pn, 
will vary. 

This is illustrated by constructing the system through 
a random deposition process, i.e. differently from the 
above procedure. The algorithm used to construct the 
system is the following. At each time step, an empty 
site and a particle orientation are chosen at random. If 
the particle can fit on this site (without overlap with 
other particles), then the site is occupied, otherwise a 
new random trial is made. This is similar to the "random 
sequential" problem often studied in the literature [[l6[ , 
here adapted to the Tetris model. 

This procedure leads to a maximum density of particles 
around Pmax ~ 0.75 above which it becomes impossible 
to add new particles. 

Differently from the previous case, in the random se- 
quential deposition simulations could not have been per- 
formed using the transfer matrix algorithm and thus we 
generated systems of size up to 500 x 1500. We studied 
the dilatancy angle in such systems stopping the con- 
struction at different p values, averaging for each p over 



100 realizations. Fig.(H) shows the estimated dilatancy 
angle which is definitely different from the data of Fig. (|^). 
In particular the onset of dilatancy is determined to be 



Pfl = 0.70 ±0.01 



(4) 
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FIG. 8. 

However, this procedure is not expected to induce long 
range correlations in the particle density or orientation, 
and thus, we expect that the universality class of the 
model remains unchanged. In particular, the critical be- 
havior Eq. His expected to hold with the same exponents. 
Although the system sizes are much smaller in the present 
case, our data are consistent with such a law. 



V. BALLISTIC DEPOSITION UNDER GRAVITY 

Finally we would like to point out another property 
related to the texture of the medium. Up to now the two 
procedures followed to generate the packing of particles 
did not single out any privileged direction. 

We now construct the packing by random deposition 
under gravity. Particles with a random orientation are 
placed at a random x position, and large y. Then the 
particle falls (along —y) down to the first site where it hits 
an overlap constraint. In this way, the packing assumes 
a well defined bulk density p fa 0.8. 

We used this construction procedure to generate lat- 
tices of size 500 x 1500 (averaged over 500 samples) cut- 
ting out the top part of the lattice which is characterized 
by a very wide interface and a non-constant density pro- 
file. On this configuration (and thus at a fixed density) 
we measured the dilatancy angle for different orientations 
of the imposed displacement on the wall with respect to 
"gravity" . 

Table 1 shows the resulting dilatancy angles obtained 
for the same density p = 0.81 using different construc- 
tions: 



• the dilution of the ordered state, 

• the sequential deposition, (in both of these cases 
the dilatancy angle does not depend on the orien- 
tation of the motion). It is worth noticing how a 
direct comparison between this case and the others 
is not possible because with the Random Sequen- 
tial Deposition one cannot obtain densities larger 
that ~ 0.75. 

• the ballistic deposition using a displacement along 
—y (against gravity), y (along gravity), and x (per- 
pendicular to gravity). In the latter case, we could 
study the problem for two orientations of the semi- 
infinite line {x = and y > or y < 0). We checked 
that the dilatancy angle was not dependent on this 
orientation. 

The data reported in Table 1 indeed shows that the 
dilatancy angle can be dependent on the direction of the 
imposed displacement. This measurement is thus sensi- 
tive to texture effects. As a side result, we note that the 
usual characterization of the dilatancy in terms of a single 
scalar (angle), albeit useful, is generally an oversimplifi- 
cation for textured media. Indeed, a number of studies 
have revealed |17[ that granular media (even consisting 
of perfect spheres) easily develop a non isotropic texture 
as can be shown by studying the distribution of contact 
normal orientations. This remark is almost obvious from 
a theoretical point of view, however, few attempts have 
been made to incorporate these texture effects in the di- 
latancy angle or even more generally in the rheology of 
granular media. 



Method 


Orientation 


iIj 


Dilution 


±x±y 


31.0 ±0.1 


BDG 


-y 


6.6 ± 0.5 


BDG 


±x 


23.8 ± 0.5 


BDG 


+y 


5.8 ± 0.5 



TABLE I. Results for the dilatancy angles obtained using 
differently prepared samples and different displacement ori- 
entations at p ~ 0.8. Dilution indicates samples obtained by 
diluting a perfect monocrystal (see text) to the desired den- 
sity; BDG indicates samples obtained with a Ballistic Depo- 
sition procedure under Gravity. We cannot compare directly 
in this table the results obtained with the Random sequen- 
tial deposition procedure (RSD) because, as mentioned in the 
text, this procedure leads to a maximum density of particles 
around Pmax ~ 0.75 above which it becomes impossible to 
add new particles. 



VI. CONCLUSION 

We have shown that dilatancy can be precisely defined 
in the Tetris model, and that it is a function of the 
density as it is well known for granular media. The onset 
of dilatancy, i.e. the "critical state" of soil mechanics, 
has been shown to corresponds to a directed percolation 
threshold density, hence justifying the term "critical" in 
this expression. Form this point of view it is important 
to stress how any comparison of our approach with real 
granular materials should be done in the neighborhood of 
the critical point where we expect a largely universal (in 
the sense of critical phenomena in the statistical physics 
vocabulary) behavior. Using different lattices we expect, 
for instance, to recover the same critical behavior (same 
exponents) but not the same values for the critical den- 
sity. To our knowledge this is the first time that such a 
mapping is proposed. We have also shown that the dila- 
tancy angle was not only determined by the density but 
also by the packing history. Finally, we have shown from 
a simple anisotropic construction that texture affects the 
dilatancy angle, even for a fixed density. 
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VII. FIGURE CAPTIONS 

Fig.(l) Schematic view of shearing of granular media in a 
shear cell. The upper part of the cell moves only 
if the medium dilates so that the direction of the 
motion forms an angle ■0, the dilatancy angle, with 
the horizontal direction. 

Fig. (2) Illustration of the Tetris model. The sites of a 
square lattice can host elongated particles shown as 
rectangles. The width and length of the particles 
induce geometrical frustrations. 

Fig. (3) Procedure used to define the dilation angle. All 
particles located on a semi-infinite line (the parti- 
cles enclosed in the round-edge rectangle on the 
left-hand side of the lattice) are moved by one 
lattice unit in the horizontal direction (shown by 
the arrows). Using the hard-core repulsion be- 
tween particles, we determine the particles which 
are pushed (shown in black) and those which may 
stay in place (grey). For each column we consider 
the lowest (in general the most external) black site 
(The gray particles within the cluster of black par- 
ticles do not play any role in the determination of 
the dilation angle). The curve connecting all these 
points defines the profile of the pushed region. The 
line connecting the first and the last points of this 
profile determines the angle, ijj, with respect to the 
direction of motion. This angle provides the value 
of the dilation angle for the particular realization 
considered. The dilation angle is actually measured 
performing an average over a large number of real- 
izations. 

Fig. (4) Lattice over which directed site percolation is tak- 
ing place. The arc bonds connect second neighbors 
along the x axis (horizontal). 

Fig. (5) Schematic representation of the mobilized region 
in the shearing procedure. Starting from p = 1, 
where one has a dilation with an angle of 7r/4, the 
dilation angle reduces until (for p — pn). A fur- 
ther reduction of p brings the system in a subcrit- 
ical regime where only a finite layer of thickness 
^// (X {pn — p)""// along the y-axis is mobilized and 
the system compactifies. 

Fig. (6) Shapes of the boundaries of two clusters for (a) 
p = 0.58 (close to the threshold for a vanishing 
dilatancy) and (b) p — 0.7. The clusters mobilized 
are above and in both cases the line interpolating 
linearly between the first and the last point the 
boundaries defines the dilation angle. 

Fig. (7) Dilatancy angle as a function of the density p in the 
case of a random dilution of the perfectly ordered 
Tetris model. The dashed line represents a fit 



obtained using Eq. (|) with pR = 0.583 ± 0.001. 
The relative errors diverge at the transition. 

Fig. (8) Dilatancy angle as a function of the density p in 
the case of a random sequential deposition. The 
dashed line represents a fit obtained using Eq. (g) 
with pn = 0.70 ± 0.01. The relative errors diverge 
at the transition. 



